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ERROR  ANALYSIS  IN  THE  JOINT  EVENT  LOCATION/SEISMIC  CALIBRATION 

INVERSE  PROBLEM 

William  Radi 
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Contract  No.  F19628-03-C-0I09 


ABSTRACT 


This  project  is  developing  new  mathematical  and  computational  techniques  for  analyzing  the  uncertainty'  in  seismic 
event  locations,  as  induced  by  observational  errors  and  errors  in  travel- time  models:  The  analysis  is  being  done  in 
the  context  of  the  multiple-event  inverse  problem,  in  which  the  locations  of  multiple  events  are  inferred  joindy  with 
travel-time  corrections  for  the  event-station  paths.  The  premise  is  that  one  of  the  events  is  the  target  of  the 
uncertainty  analysts  while  the  others  are  calibration  events  that  serve  to  constrain  the  path  corrections  within  some 
level  of  error.  An  uncertainty  analysis  on  the  coupled  location/calibration  inverse  problem  leads  to  the  notion  of  a 
“multiple-event”  confidence  region  on  the  target  event  location,  which  accounts  implicitly  for  errors  in  the  inferred 
path  corrections,  including  the  effects  of  uncertainty  in  the  calibration  event  locations.  The  approach  we  have 
developed  considers  Ihe  nonlinearity  of  the  forward  problem  and  allows  the  use  of  non -Gaussian  models  for  data 
errors  and  soft  or  hard  constraints  on  the  problem  parameters. 

The  project  has  addressed  mainly  the  basic  multiple- event  location  problem,  wherein  travel- time  corrections 
comprise  a  simple  time  term  for  each  station/phase  combination  in  the  data  set.  We  have  developed  a  numerical 
scheme  for  computing  multiple-event  confidence  regions  for  this  problem,  based  on  grid  search  and  Monte  Carlo 
sampling  techniques.  Given  the  computationally  intensive  nature  of  the  approach,  some  approximations  were 
developed  to  reduce  the  computational  effort.  Experiments  with  data  from  the  Nevada  Test  Site  validated  the 
approach  and  demonstrated  the  adequacy  of  the  approximations  when  the  location  uncertainty  is  not  too  large.  A 
modification  of  the  approach  we  are  investigating,  which  offers  the  possibility  of  additional  efficiency,  is  the 
computation  of  multiple-event  confidence  regions  defined  in  a  Bay  esian  sense.  Moreover,  we  are  investigating  an 
extension  of  our  multiple-event  location  algorithm  to  the  case  where  travel-time  corrections  are  parametrized  as  a 
function  of  the  event  location 
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OBJECTIVE 


This  project  is  developing  new  mathematical  and  computational  techniques  lor  quantifying  the  errors  in  seismic 
event  locations,  including  the  effects  of  observational  errors  and  errors  in  the  travel-time  forward  model  Our 
approach  associates  the  latter,  or  model  errors,  with  the  uncertainty  in  path  travel- time  corrections  that  have  been 
inferred  from  calibration  data  We  thus  analyze  event  location  uncertainty  in  the  context  of  the  joint 
location/calibration  inverse  problem;  using  arrival  time  data  from  multiple  events  and  stations  to  simultaneously 
locate  the  events  and  estimate  calibration  parameters  that  determine  the  travel- time  corrections.  Calibration 
parameters  can  be  as  simple  as  lime  terms  at  stations  or  as  complex  as  a  3D  velocity'  model.  Our  approach  treats  one 
of  the  multiple  events  as  a  new  event  under  investigation  and  the  remaining  events  as  calibration  events  Performing 
uncertainty  analysis  on  the  joint  location/calibration  problem  accounts  for  the  key  sources  of  error  in  the  new  event 
location,  including  picking  errors  in  the  arrival  times  for  the  new  and  calibration  events  and  errors  in  the  calibration 
event  locations.  We  are  addressing  this  joint  inverse  problem  with  numerical  techniques,  like  grid  search  and 
Monte-Carlo  simulation,  that  avoid  the  limitations  of  analytic  approaches  to  uncertainty  analysis  in  large  inverse 
problems.  These  limitations  include  the  restriction  to  Gaussian  data  errors,  the  necessity  of  using  a  linear 
approximation  to  the  travel-time  forward  model,  and  restricted  mechanisms  for  incorporating  a  priori  constraints  on 
the  unknowns. 

In  previous  years'  efforts,  we  implemented  and  tested  this  approach  to  event  location  uncertainly  for  a  relatively 
simple,  but  widely  used,  parameterization  tor  travel-time  corrections  (time  terms).  These  efforts  demonstrated  the 
conceptual  validity  of  the  approach,  w'hich  could  be  validated  quantitatively  in  simple  situations  but  which  revealed 
the  practical  difficulties  associated  with  its  computational  intensity.  Here  we  report  on  some  incremental 
improvements  to  the  uncertainty  algorithm  we  have  developed. 

RESEARCH  ACCOMPLISHED 


Summary  of  Approach 

Rodi  (2004,  2005)  presented  the  mathematical  formulation  of  the  joint  location/calibration  inverse  problem  and 
maxi  mum -likelihood  uncertainty'  approach  that  are  the  basis  for  this  project.  To  summarize,  we  write  the  joint 
inverse  problem  as 


dij  =  Ti(xi)  +  ti+cij+^,  (1) 

where  i  indexes  each  of  m  seismic  events  and  j  indexes  each  of  n  station/phase  combinations  that  have  been 
observed  from  one  or  more  of  the  events.  Then,  du  denotes  the  arrival-time  observation  lor  the  j'th  event  and >th 
station/phase  ((i/)th  path),  x,  and  tf  are  the  origin  parameters  (hypocenter  and  time,  respectively)  of  the  ith  event,  7 ) 
is  a  model-based  travel-time  function  for  theytb  station/phase,  ctj  is  a  correction  to  this  function  for  path  <7 J)f  and  etJ 
is  an  observational  error.  This  equation  holds  only  for  the  paths  (ij)  for  which  data  have  been  observed  The  joint 
location/calibration  problem  is  often  referred  to  as  the  multiple-event  location  problem  The  project  has  focused  on 
the  4  basic  multiple-event  location  problem  w'hereby  the  path  corrections  are  assumed  to  be  event- independent  (see 
Jordan  and  Sverdrup,  1 98 1 ,  and  Pavlis  and  Booker,  1 983): 

S  =  V  (2) 

The  unknown  parameters  of  the  inverse  problem  are  then  the  event  hvpocenters  and  origin  times,  %h  /„/=!,.  jti, 
and  the  time  terms  for  the  stalion/phase  combinations;  ajJ  =  1T,.. ,  n 

The  premise  of  our  approach  is  that  one  of  the  events  (e.g.,  i=  1)  is  the  target  of  the  location  uncertainly  analysis, 
w'hile  the  remaining  events  (i  =  2>...  jn)  are  calibration  events,  whose  data  provide  information  about  the  ar 
Generally,  prior  information  is  available  to  constrain  the  location  of  one  or  more  of  the  calibration  events 

The  likelihood  function  for  the  multiple-event  location  problem  is  determined,  in  part,  by  the  assumed  probability 
distribution  ot  the  observational  errors  in  the  data,  etj.  We  take  this  to  be  a  generalized  Gat4ssiatt  distribution  of  order 
p  (Billings  cl  a  I,  1 994).  Additionally,  the  likelihood  function  depends  on  the  prior  information  available  for  the 
calibration  event  locations.  Rodi  (2005)  discusses  the  use  of  hard  and  soft  prior  constraints  on  an  event  location 
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The  latter,  which  we  consider  in  this  paper,  take  the  form  of  a  prior  probability  distribution  on  the  event  hypocenter. 
As  with  the  data  errors,  we  allow  this  distribution  to  be  of  the  generalized  Gaussian  type. 

Given  these  assumptions,  the  likelihood  1  unction  for  the  basic  multiple-event  location  problem  is  given  by 
-logL  =  conrt+lj-Ul,  +2  1  |K  -  l0i|T 

where  ||  ||  denotes  the  Euclidean  distance  in  space.  The  prior  information  on  the  ith  event  location  is  specified  by  a 
prior  location  estimate,  x0/;  a  standard  error,  oQi ;  and  a  generalized  Gaussian  order,  p0h  where  we  generally  assume 
Poj  >  2.  The  generalized  Gaussian  order  for  the  data  errors  is p(l  <p<2\  and  the  data  standard  errors  are  aQi.  For 
the  purpose  of  exposition,  we  assume  the  data  standard  errors  are  known 

Location  Confidence  Regions 

Our  approach  to  uncertainty  analysis  m  the  joint  location/calibration  inverse  problem  separates  the  unknown 
parameters  into  two  vectors:  p,  containing  the  subset  of  “target”  parameters  on  which  a  confidence  region  is 
desired— and  q,  containing  the  remaining  “nuisance'*  parameters.  For  example,  if  p  -  (the  hypocenter  of  the 
target  event)  then 


(tli  Ijj  fjt-  -  ,  Ini  iflB)  -  (4) 

To  consider  a  confidence  region  on  the  epicenter  of  the  target  event,  we  would  move  its  depth  parameter,  zu  from  p 
to  q.  Given  this  separation  of  parameters,  we  will  show  the  likelihood  function  as  L(p,q;d),  where  d  denotes  the 
vector  of  arrival* time  observations,  dtJ. 

A  confidence  region  on  p  is  determined  by  the  behavior  of  the  likelihood  function  L  as  a  function  of  p, 
encompassing  the  points  that  yield  the  largest  values  of  L.  Throughout  this  project,  vve  have  defined  such  regions 
using  the  approach  of  Neyman- Pearson  (N-P)  hypothesis  testing.  We  take  as  a  test  statistic,  T(p,d>,  the  log  likelihood 
ratio  given  by 


r(p,d)  =  log  max  f(pTq;d)  -  log  max  L(p.qd) 

M  q  (5) 

■  log  L(p\q  *;d)  -  log  max  £.(p,q;d). 

q 

where  (p*7q*>  denotes  the  maximum-likelihood  solution  for  the  parameters.  We  see  that  rcompares  the  likelihoods 
that  are  achieved  with  p  fixed  to  a  particular  value  (second  term  of  t)  and  with  p  free  to  vary  (first  term).  In  each 
case  q  is  free.  A  N-P  confidence  region  on  p  at  confidence  level  ft  (e.g.  ft  =  95%)  is  the  locus  of  points  satisfying 

*(ivd)<r^ph  (6) 


where  r^p)  is  a  critical  value  of  the  probability'  distribution  of  r,  as  induced  by  the  distribution  of  the  data  errors: 

Prob  {z^p,  d)  <  r^p)  |  p}  =  ft.  (7) 

In  writing  this  and  Equation  (6),  we  have  allowed  for  the  possibility'  that  the  distribution  of  r  may  depend  on  p,  but 
vve  have  tacitly  assumed  that  is  does  not  depend  on  the  nuisance  parameters. 

We  have  developed  a  two-step  numerical  algorithm  to  implement  the  N-P  approach  to  confidence  regions.  The  first 
step  maps  the  test-statistic,  i(p,  d),  on  a  grid  in  p-space.  This  entails  minimizing  the  likelihood  function  with  respect 
to  q  with  p  fixed  in  turn  to  each  point  on  the  grid.  The  second  step  performs  Monte  Carlo  simulation  to  inter  the 
critical  values  r^p).  This  entails  maximizing  L  with  respect  to  q,  with  p  fixed  and  free,  for  each  of  many  realizations 
ot  synthetic  data.  Both  steps  of  this  algorithm  are  computationally  intensive  because  they  involve  repeated 
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maximization  ot  L,  i.e.,  repeated  solution  ol  a  multiple-event  location  problem  involving  many  location  and 
calibration  parameters. 

Needless  to  say*  this  formulation  applies  to  the  single-event  location  problem  and  "  single-event  confidence  regions" 
by  removing  the  travel- time  corrections,  aJy  and  the  calibration  events  and  their  data  from  the  problem 
(i.e.,a,  =  0,m  =  1). 


Monte  Carlo  Short  Cuts 

When  data  and  prior  enors  are  Gaussian  {p  =p0i  =  2)  and  when  the  model-based  travel-time  functions  can  be 
approximated  as  linear,  t q  does  not  depend  on  p,  In  general,  however,  Tfi  might  depend  on  p,  indicating  the  need  to 
pertorm  a  separate  Monte  Carlo  simulation  with  each  value  of  p  on  the  sampling  grid  used  for  the  “true"  parameters 
in  calculating  synthetic  data  realizations.  This  is  not  practical  for  realistic  problems,  and  our  algorithm  performs 
Monte  Carlo  simulation  only  with  p  set  to  its  maximum-likelihood  estimate ,  p*  That  is.  Equation  (6)  is  replaced 
with 


(8) 

This  approximation  has  been  assumed  by  others  who  have  computed  numerical  confidence  regions  on  seismic  event 
locations  (e.g.,  Wilcock  and  Toomey,  1991 ).  However,  to  the  extent  that  the  dependence  on  p  exists,  the  confidence 
regions  based  on  Equation  (8)  lose  the  "confidence  property,”  i.e.,  the  region  includes  the  true  parameter  vector  with 
probability  regardless  of  what  the  true  parameter  vector  is. 

A  second  significant  short  cut  we  have  implemented  pertains  to  the  nuisance  parameters.  In  the  Gaussian/] inear  case, 
the  probability  distribution  of  r<  which  is  chi-squared)  depends  only  on  the  dimensionality  of  p,  regardless  of  what, 
if  any,  nuisance  parameters  there  are.  As  a  result,  q  can  be  ignored  for  the  purpose  of  determining  t In  a  Monte 
Carlo  simulation  of  this  implies  that  the  likelihood  maximization  with  synthetic  data  realizations  can  treat  q  as 
fixed.  In  the  location/calibration  problem,  most  of  the  computational  benefit  of  doing  this  is  achieved  by  fixing  the 
calibration  event  locations  so  that  each  likelihood  maximization  involves  the  work  of  single-event,  rather 
than  multiple-event,  location 

Quasi- Bayesian  Confidence  Regions 


In  Bayesian  inference,  the  likelihood  function  we  defined  m  Equation  (3),  including  the  factor  for  prior  information 
on  the  parameters,  is  treated  as  an  unnormalized  posterior  probability  density'  function  on  the  parameters.  A 
Bayesian  confidence  region  on  p  at  confidence  level  fS  is  a  region  in  p- space,  P^h  which  includes  a  proportion  of 
the  posterior  probability  density; 

,  ({P Jyi  rfq/-(p,q.d)  =  ^//,|  dqL( p,q,d) ,  (9) 


where  Pr  and  Q ,  denote  the  entire  p-space  and  q-spaee,  respectively  In  a  Gaussian/linear  problem,  Bayesian 
confidence  regions  coincide  with  N-P  regions  and  thus  posses  the  confidence  property.  In  genera],  however,  the  two 
types  of  regions  differ,  and  Bayesian  ones  do  not  necessarily  have  the  confidence  property. 

1  he  difference  between  N-P  and  Bayesian  confidence  regions  in  the  presence  of  nonlinearity  is  illustrated  in 
Figure  1  This  is  a  “toy"  problem  involving  one  unknown  parameter — event  depth — and  a  single  datum  that  is  taken 
to  be  a  direct  estimate  of  the  depth.  We  introduce  nonlinearity  into  the  problem  in  the  form  of  the  hard  constraint 
that  the  true  event  depth  (horizontal  axis)  cannot  be  negative,  even  though  the  depth  estimate  (vertical  axis)  might 
be,  such  as  when  arrival  times  yield  an  "airquake"  as  the  best-fitting  solution  for  an  event  location.  Intuitively,  the 
positivity  constraint  on  depth  should  affect  a  confidence  region  on  depth  when  the  depth  estimate  is  near  z  =  0T  and 
we  can  see  from  the  figure  that  it  does.  The  case  in  the  left  panel  assumes  that  the  depth  estimate  has  a  Gaussian 
distribution  [p-  2)  while  the  case  on  the  right  assumes  a  Laplace  distribution  {p  -  1  y  The  green  lines  in  each  plot 
trace  the  confidence  limits  for  a  N-P  confidence  intervals,  which  possess  the  confidence  property.  The  red  dots  are 
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the  approximate  N-P  confidence  intervals  that  result  when  the  critical  value  Tpis  considered  only  for  the  estimated 
(maximum-likelihood)  value  of  depth,  as  per  Equation  (8),  and  the  blue  dots  mark  the  Bayesian  confidence  limits. 
Generally,  the  approximate  N-P  intervals  are  smaller  than  the  Bayesian  intervals  for  shallow  depth  estimates, 
showing  they  are  more  influenced  by  the  positivity  constraint.  Which  one  better  represents  the  exact  N-P  intervals 
depends  on  the  type  of  data  distribution  used. 

True  Depth  (km)  True  Depth  (km) 


Figure  L  Confidence  intervals  on  event  depth  in  a  fictitious  problem  involving  a  direct  estimate  of  depth, 

which  is  either  Gaussian  (left  panel)  or  Laplace  (right)  distributed,  with  a  standard  error  of  10  km  in 
each  case.  Three  types  of  confidence  limits  (at  0  =  90%)  are  show  n,  each  as  a  function  of  the  depth 
estimate  itself  (vertical  axis):  exact  N-P  (green  line),  approximate  N-P  (red  dots),  and  Bayesian  (blue 
dots).  The  corresponding  confidence  interv  al  is  the  range  of  true  depth  between  these  limits. 


It  is  straightforward  to  modify  our  numerical  algorithm  to  compute  Bayesian  confidence  regions  if  we  make  an 
additional  approximation,  namely,  replacing  the  integration  over  the  nuisance  parameters  q  in  Equation  (9)  with 
maximization.  (Acknowledging  this  approximation,  we  will  call  the  resulting  confidence  regions  “quasi- Bayesian”) 
Considering  the  definition  of  r  in  Equation  (5),  Equation  (9)  becomes  equivalent  to  Equation  (6),  with  a  single 
recomputed  as  the  solution  of 
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Application  to  NTS  Data 

We  are  testing  our  new  uncertainty  analysis  techniques  on  regional  seismic  arrival  times  from  NTS  explosions,  for 
which  precise  locations  and  origin  times  are  known  (Walter  et  at  ,  2003),  The  results  shown  here  are  based  on  a 
subset  of  the  available  data  set,  comprising  262  Pn  arrivals  from  33  explosions  at  Fahute  Mesa  and  Rainier  Mesa. 
This  subset  included  arrivals  from  the  27  stations  that  spanned  the  epicentre!  distance  range  l,5°-9°  and  had 
Pn  arrivals  for  at  least  2  of  the  33  events. 

When  computing  a  confidence  region  on  the  location  of  one  of  the  events  (the  “target”  event),  the  remaining 
32  events  were  treated  as  calibration  events.  Only  one  calibration  event  was  assigned  a  finite  ground  truth  (GT) 
level;  a  relatively  we  11- recorded  Rainier  Mesa  explosion  (16  Pn  arrivals).  Its  assigned  GT  accuracy  w  as  varied 
between  0  and  5  km.  The  depths  of  all  events  were  fixed  to  their  true  values  in  these  tests.  The  origin  times  of  all 
events,  including  the  GT  calibration  event,  were  unconstrained.  The  IASP91  travel-time  tables  were  used  for  die 
forward  model. 
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figure  2  compares  the  epicenter  solutions  obtained  with  single-event  location  (zero  travel-time  corrections)  with 
multiple-event  solutions  obtained  with  varying  GT  levels  for  the  GT-cahbration  event.  We  see  that  the  multiple- 
event  solutions  achieve  greater  accuracy  in  the  relative  locations  between  events,  but  the  absolute  event  mislocations 
become  larger  as  the  location  constraint  on  the  GT- calibration  event  is  relaxed. 
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Figure  2*  Four  epicenter  solutions  for  33  NTS  explosions  in  the  Pahute  Mesa  and  Rainier  Mesa  testing  areas. 

Each  solution  location  (blue  dot)  is  connected  by  a  line  to  the  corresponding  true  (GTO)  epicenter  of 
the  event  (red  circle).  Top  left:  Solution  obtained  with  single-event  location  (travel-time  corrections 
set  to  zero).  Top  right:  Solution  obtained  with  multiple-event  location  with  the  location  of  one 
well-recorded  Rainier  Mesa  event  (37.21  °N,  1 16.21°W)  held  fixed  (GTO  constraint).  Bottom: 
Multiple-event  location  solutions  with  the  same  Rainier  Mesa  event  treated  as  GT2  (bottom  left)  and 
GT5  (bottom  right). 

Single-Event  and  Multiple-Event  Confidence  Regions 

figure  3  compares  single-event  and  multiple-event  confidence  regions  for  a  Pahute  Mesa  event  with  relatively  few 
arrivals  (6  Pn  arrivals).  The  regions  were  calculated  using  a  Gaussian  distribution  for  observ  ational  errors  {p  -  2). 
The  colored  areas  are  the  confidence  regions,  for  different  confidence  levels,  computed  with  our  two-step 
(approximate  N-P)  numerical  algorithm.  The  single-event  regions  (left)  do  not  include  the  effect  of  errors  in 
travel-time  corrections  (model  errors)  and  agree  almost  exactly  with  analytically  derived  confidence  ellipses  based 
on  only  the  picking  error  variance  (shown  for  /?  =  95%). 
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Figure  3.  Epicenter  confidence  regions  for  a  Pahute  Mesa  explosion  with  6  Pn  arrivals,  computed  with  our 
two-step  numerical  algorithm.  Left:  Single-event  confidence  regions  (path  travel-time  corrections 
assumed  known).  Right:  Multiple-event  confidence  regions,  treating  a  well-recorded  Rainier  Mesa 
event  as  a  GT0  calibration  event.  Confidence  regions  are  shown  for  90%,  95%,  and  98%  confidence 
(blue,  green,  and  red,  respectively).  In  each  frame  the  black  circle  marks  the  maximum-likelihood 
solution  for  the  event  location,  and  the  white  circle  is  its  GTO  location  (from  Walter  et  al.,  2003).  The 
error  distribution  for  picking  errors  w  as  assumed  to  be  Gaussian  (p  =  2).  The  ellipse  in  each  frame  is 
the  conventional  95%  single-event  confidence  ellipse,  computed  analytically  under  the 
Gaussian/linear  assumption. 

The  multiple-event  confidence  regions  (right  pane!)  are  larger  because  they  account  for  the  additional  location 
uncertainty  induced  by  uncertainty  in  the  estimated  travel-time  corrections.  The  multiple-event  confidence  regions 
in  this  case  assumed  that  the  G 1  event  at  Rainier  Mesa  was  GTO.  Therefore,  the  confidence  regions  do  not  account 
for  any  location  uncertainty  in  this  calibration  event. 

Monte  Carlo  Simulation  vs.  Likelihood  Integration 


Figure  4  shows  multiple-event  confidence  regions  for  the  same  Pahute  Mesa  target  event  just  considered,  computed 
two  different  ways.  The  results  in  the  top  row'  used  the  same  two-step  procedure  used  for  Figure  3  (approximate 
confidence  regions).  The  results  in  the  bottom  row  are  quasi- Bayesian  (q-B)  confidence  regions.  The  two  types 
of  confidence  regions  are  generated  from  the  same  numerical  mapping  of  the  test  statistic  (t)  as  a  function  of 
epicenter,  but  they’  use  different  critical  values  of  the  test  statistic  (rs)  to  determine  which  contour  of  t  corresponds 
to  each  confidence  level.  The  N-P  regions  determine  rf  with  Monte  Carlo  simulation  (doing  500  realizations),  while 
the  q-B  regions  use  likelihood  integration,  as  discussed  earlier. 

The  lett  panels  of  Figure  4  assume  a  GTO  calibration  event  (the  top  left  panel  repeats  the  right  panel  of  Figure  3.) 
The  center  and  right  panels,  respectively,  assume  the  GT- calibration  event  is  GT2  and  GT5.  We  see  that  the 
confidence  regions  grow  as  the  G  f  level  is  increased  (see  Rodi  (2005)  for  additional  examples).  Comparing  the  two 
methods  (top  vs.  bottom  panels),  we  see  that  the  quasi-Bayesian  confidence  regions  are  smaller.  This  owes  to  the 
fact  that  likelihood  integration  yielded  smaller  values  of  Te  than  did  Monte  Carlo  simulation. 

T  able  1  compares  the  critical  statistics  determined  by  the  two  methods  with  the  value  predicted  by  the 
Gaussian/linear  theory.  The  critical  values  from  likelihood  integration  and  the  linear/Gaussian  theory  agree  in  both 
the  single-event  and  multiple-event  cases,  but  the  Monte  Carlo  values  are  larger  than  either  in  the  multiple -event 
cases.  A  discrepancy  from  the  Gaussian/linear  values  can  only  be  due  to  nonlinearity  (since  Gaussian  data  and  prior 
errors  were  assumed).  However,  we  w'ould  expect  the  effect  of  nonlinearity  to  increase  as  the  GT  level  of  the 
GT -calibration  event  increases,  but  the  Monte-Carlo  values  do  not  increase.  Based  on  this  comparison,  we 
tentatively  prefer  the  likelihood-integration  method  for  determining  The  remaining  examples  of  this  paper  show 
quasi-Bayesian  confidence  regions. 
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Table  1.  Comparison  of  critical  %  values  from  Monte  Carlo  simulation  (M.C),  likelihood  integration  (L.L) 
and  Gaussian/li near  theory  \  j 
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Figure  4.  Approximate  N-P  (top)  and  quasi- Bayesian  (bottom)  multiple-event  confidence  regions  for  the  same 
Pahute  Mesa  explosion  considered  in  Figure  3,  The  Rainier  Mesa  GT-ea  lib  ration  event  is  assumed  to 
be  GTG  (left),  GT2  (center),  or  GTS  (right).  A  Gaussian  error  model  is  used  for  both  the  picking 
errors  and  the  prior  location  error  ip  =/?0  =  2).  Plotting  conventions  are  the  same  as  Figure  3* 

Non-Gaussian  Error  Models 

The  examples  shown  above  assumed  Gaussian  data  errors  (p  =  2)  and  a  Gaussian  prior  error  on  the  GT-cahbration 
event  at  Rainier  Mesa  (p0  =  2).  Using  other  tests  with  NTS  data,  Rodi  (2005)  showed  that  increasing  p0  ("hardening^ 
the  GT  constraint)  caused  confidence  regions  to  become  non-elliptical,  especially  in  the  GTS  case.  Figure  5 
illustrates  this  for  the  target  event  we  have  been  considering, 

Our  final  test  shows  the  effect  of  using  a  non-Gaussian  probability  distribution  for  the  observational  errors  Figure  6 
shows  single-event  and  multiple- event  confidence  regions  (same  target  event)  for  two  values  of  the  generalized 
Gaussian  distribution;  p-  1.5  (top)  and  1  25  (bottom).  The  confidence  regions  depart  from  an  elliptical  shape  more 
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lorp  1 .25  (bottom  panels)  than  for p  —  1.5,  as  expected.  The  non-ellipticity  is  most  noticeable  in  the  single-event 
regions,  which  are  most  strongly  influenced  by  the  picking  error  distribution. 

These  results  also  show  the  numerical  challenge  of  automatically,  but  efficiently,  grid  ding  the  test  statistic  when  its 
behavior  departs  significantly  from  a  quadratic  function  (see  lower  panels).  Our  current  scheme  designs  a  mapping 
grid  in  advance  (guided  by  the  conventional  error  ellipse  parameters),  but  these  tests  indicate  the  need  for  an 
adaptive  gridding  scheme. 


calibration  event:  generalized  Gaussian  orders  p„  =  2,  5,  and  20  (left,  center,  and  right,  respectively). 
The  GTS  constraint  is  at  90%  confidence  in  each  case.  The  data  errors  are  assumed  to  he  Gaussian. 
The  target  event  is  the  same  as  in  previous  figures. 
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Figure  6:  Single-event  (left)  and  multiple-event  confidence  regions  computed  using  two  non-Gaussian 
probability'  distributions  for  picking  errors:  generalized  Gaussian  orders  p  =1.5  (top)  and 
1.25  (bottom).  The  prior  error  on  the  GT-ca  lib  ration  event  for  the  multiple-event  regions  is  either 
GTO  (center  panels)  or  GT2  (right),  with  a  Gaussian  prior  error  assumed  for  the  latter. 
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CONCLUSIONS  AND  RECOMMENDATIONS 


The  project  has  made  significant  progress  in  the  development  of  a  new  approach  to  event  location  uncertainty 
analysis.  The  approach  is  based  on  the  framework  of  the  joint  locatkm/calibration  inverse  problem,  which 
considers  the  effects  of  both  picking  and  model  errors  on  event  location  errors.  We  have  implemented  numerical 
techniques  for  computing  "multiple-event”  confidence  regions  without  need  of  the  conventional  assumptions  of 
Gaussian  errors  and  linear ity  of  the  forward  problem.  Our  recent  work  has  focused  on  additional  enhancements  to 
the  computational  efficiency  of  the  approach,  which  have  led  to  the  alternative  of  using  likelihood  integration  to 
replace  the  Monte  Carlo  simulation  step  of  our  original  algorithm.  This  alternative  is  tantamount  to  a  Bayesian 
approach  to  uncertainty  but  can  also  be  viewed  as  just  a  different  computational  short  cut  for  approximating  the 
exact  confidence  regions  defined  by  N-P  hypothesis  testing.  Our  results  to  dale  suggest  that  this  new 
approximation  may  in  fact  be  more  accurate.  It  is  apparent,  however,  that  further  work  is  needed  to  produce 
reliable  calculations  when  error  distributions  depart  significantly  from  Gaussian  and  confidence  regions  from 
ellipses. 

The  practicality  of  our  approach  for  more  comptex  parameterizations  of  travel-time  corrections,  such  as  velocity 
models,  is  in  question  at  this  point.  However,  our  future  plans  include  pursuing  a  more-efficient  scheme  for 
generating  likelihood  functions  on  a  grid  (e.g.,  using  linearization  w  ith  respect  to  "nuisance"  parameters),  and 
this  may  alter  the  outlook  In  any  event,  extending  our  approach  to  more-complex  parameterizations  remains  a 
feasible  task  for  research  that  we  intend  to  pursue  during  the  final  year  of  the  project. 
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